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A  NONLINEAR  HIERARCHICAL  MODELLING  APPROACH  FOR  CENSUS 

UNDERCOVERAGE  ESTIMATION 

Yong  You1 

ABSTRACT 

Area-level  nonlinear  mixed  effects  models  are  considered  in  this  paper  for  Canada  census  undercoyerage  estimation. 
We  fit  an  area-level  nonlinear  mixed  effects  model  to  the  province-level  undercoverage  survey  estimates.  In  particular, 
the  sampling  model  is  based  on  the  survey  estimate  of  the  undercoverage  count,  and  the  linking  model  is  a  log-linear 
model  for  the  undercoverage  rate.  A  full  hierarchical  Bayes  (HB)  approach  is  developed  to  obtain  the  posterior 
estimates  of  the  census  undercoverage  using  Markov  Chain  Monte  Carlo  (MCMC)  sampling  methods.  Our  result 
shows  that  the  proposed  method  can  provide  efficient  model-based  estimates.  Analysis  of  model  fitting  is  also 
presented  using  posterior  predictive  distributions,  and  the  corresponding  result  indicates  that  the  proposed  model  fits 
the  data  quite  well. 

KEY  WORDS:  Census  undercoverage,  Gibbs  sampling,  Model  checking,  Nonlinear  mixed  model,  Posterior. 

RESUME 

Les  modeles  non  linSaires  a  effets  mixtes  au  niveau  des  regions  sont  considered  dans  cet  article  dans  le  cadre  de 
l'estimation  de  la  sous-couverture  pr6sente  au  recensement  canadien.  Un  modele  non  lineaire  a  effets  mixtes  au  niveau 
des  r6gions  est  utilise'  afin  de  moddliser  les  estimes  de  sous-couverture  au  niveau  provincial.  En  particulier,  le  modele 
d'echantillonnage  est  bas6  sur  une  estimation  de  la  sous-couverture,  et  le  modele  de  liaison,  servant  a  modeliser  le  taux 
de  sous-couverture,  est  un  modele  log-lin6aire.  Une  approche  hierarchique  de  Bayes  est  developpee  afin  d'obtenir  des 
estimateurs  a  posteriori  de  la  sous-couverture,  en  utilisant  les  m&hodes  de  Monte  Carlo  avec  chaine  de  Markov 
(MCCM).  Les  r^sultats  obtenus  indiquent  que  la  mSthode  proposee  engendre  des  estimateurs,  bas6s  sur  le  modele,  qui 
sont  efficaces.  De  plus,  une  evaluation  du  modele,  utilisant  les  distributions  predictives  a  posteriori,  est  presentee  et  les 
rgsultats  indiquent  que  le  modele  propose  semble  adequat. 

MOTS  CLES:  sous-cou venire  au  niveau  du  recensement,  echantillonneur  de  Gibbs,  evaluation  du  modele,  modele  non 
lineaire  a  effets  mixtes,  a  posteriori. 


1.  INTRODUCTION 

In  Canada,  a  census  is  conducted  every  five  years. 
However,  the  census  does  not  enumerate  all  the 
inhabitants  that  should  fill  a  census  form  on  Census 
Day.  In  the  1991  Canadian  census,  it  is  estimated  that 
about  3%  of  the  population  were  not  enumerated. 
Thus  the  census  needs  to  be  adjusted  for 
undercoverage  in  order  to  properly  represent  the 
demographic  picture  of  the  country  on  Census  Day. 
Since  1966,  the  Reverse  Record  Check  (RRC)  has 
been  used  by  Statistics  Canada  to  measure  the  gross 
number  of  persons  missed  by  the  census.  Starting 
1991,  an  Overeoverage  Study  was  conducted  to 
measure,  the  gross  number  of  persons  erroneously 
included  in  the  census.  In  1991,  for  the  first  time,  the 
RRC  results  were  combined  with  those  of  the 
Overeoverage  Study  to  produce  the  direct  survey 
estimates  of  the  net  undercoverage  for  the  nation  and 
all  provinces.    Through  the  analysis  of  the  results  of 


these  coverage  studies,  the  collection  methodology  is 
adjusted  in  order  to  improve  coverage  in  the 
succeeding  census.  More  details  on  the  coverage 
studies  can  be  found,  for  example,  in  Germain  and 
Mien  (1993). 

In  1991,  the  population  estimates  were  based  on  the 
census  counts  adjusted  for  the  estimated  net 
undercoverage  in  the  census.  The  base  population  was 
formed  by  adding  the  net  provincial  undercoverage 
estimate  to  the  provincial  census  count.  This  created 
an  adjusted  base  upon  which  all  the  other  population 
figures  were  v  derived  using  modelling  and 
demographic  methods.  Rivest  (1995)  proposed  a 
composite  estimator  to  estimate  the  provincial 
undercoverage  using  the  national  undercoverage  rate 
as  a  synthetic  estimate.  The  effect  of  the  composite 
estimator  is  to  shrink  all  provincial  rates  to  the  national 
rate.  Rivest's  composite  estimator  performs  poorly  at 
extreme  provinces,  namely,  P.E.I.  and  Ontario,  the 
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smallest  and  the  largest  provinces  of  Canada  (You, 
1997). 

In  recent  years,  modeling  techniques  have  been 
applied  to  improve  the  direct  undercoverage  estimates 
from  sample  surveys.  Datta  et  al.  (1992)  and 
Zaslavsky  (1993)  proposed  hierarchical  Bayes  (HB) 
models  for  the  adjustment  of  census  undercoverage 
based  on  US  census  count  and  situations.  Dick  (1995) 
proposed  a  regression  model  utilizing  an  empirical 
Bayes  (EB)  method  to  estimate  the  undercoverage  in 
small  domains.  Dick  and  You  (1997)  employed  a 
linear  Fay-Herriot  model  (Fay  and  Herriot,  1979)  for 
census  undercoverage  rates  in  a  HB  framework.  In 
this  paper  we  are  particularly  interested  in  modeling 
the  province-level  census  undercoverage  counts.  In 
particular,  we  propose  a  nonlinear  mixed  model  in 
which  the  sampling  model  is  based  on  the 
undercoverage  counts  whereas  the  linking  model  is 
based  on  rates.  The  undercoverage  rate  is  expressed  in 
terms  of  undercoverage  counts.  The  proposed  model 
overcomes  some  limitations  of  the  model  proposed  by 
Dick  and  You  (1997)  and  can  provide  model-based 
estimates  for  undercoverage  counts  and  rates 
simultaneously.  The  remainder  of  the  article  is  laid 
out  as  follows.  In  Section  2  we  formally  present  the 
models  for  census  undercoverage  estimation.  In 
Section  3  we  illustrate  the  Gibbs  sampling  method 
based  on  the  proposed  nonlinear  hierarchical  models. 
In  Section  4  we  describe  the  implementation  of  the 
Gibbs  sampler  and  present  the  estimation  results.  In 
Section  5  we  describe  diagnostics  used  to  assess  model 
fitting.  Finally,  we  present  some  concluding  remarks 
in  Section  6. 

2.  UNDERCOVERAGE  MODELS 

Suppose  there  are  m  provinces.  In  the  i-th  province 
the  census  has  counted  c,- persons  while  an  unknown 
number  u,  persons  are  missed  by  the  census.  The 
coverage  studies  provide  a  survey  estimate  y,-  of  the 
net  undercoverage  along  with  an  associated  variance 
£f .  The  true  population  of  the  i-th  province  can  be 
written  as  Tt  =  c,.,+  u,.  Since  the  census  count  c,is 
observed  without  sampling  error,  most  of  the  work  in 
constructing  a  model  centres  on  the  estimate  of  missed 
persons  m,  .  The  coverage  studies  use  sample  surveys 
which  through  standard  estimation  procedures  produce 
an  estimate  of  the  number  of  missed  persons.  The 
model  generally  used  to  describe  this  estimation 
situation  can  be  written  as         > 

y,.  =  u,  +  £■,   h-i,...;m,        "■■  (1) 

where'' e,  -  N(0,^) .  This' model  assumes  that  the 
direct  survey  estimator  y(  is  design-unbiased.     The 


design-unbiased  assumption  may  be  restrictive,  in 
particular  the  estimates  of  missed  persons  y,  are 
probably  subject  to  some  unknown  bias  (Zaslavsky, 
1993).  However,  since  we  don't  have  estimate  for  the 
possible  bias  as  in  the  US  case  (Zaslavsky,  1993),  we 
assume  that  the  survey  estimate  y,  is  design-unbiased 
for  m,  .  This  is  also  a  common  assumption  used  in 
small  area  estimation  (Rao,  1999).  The  normal 
distribution  assumption  of  the  sampling  errors  at  the 
province  level  seems  quite  reasonable'"due  to  large 
sample  size.  The  sampling  variances  £2are  estimated 
through  generalized  variance  function  models  of  the 
form  Vr(y,)ocC|?'  (Dick,  1995)  and  then  treated  as 
known  in  the  model.  Thus  £2is  actually  a  smoothed 
estimate  of  the  sampling  variance  of  y, . 

Instead  of  using  (1)  directly  to  model  the 
undercoverage  count  u, ,  Dick  and  You  (1997)  used  a 
transformation  and  considered  modelling 
undercoverage  rate.  Namely,  let  di  =  w, /(u, ,+  c, )  and 
r.  =y. /(y.  +c,),  where  0,-is  defined  as  the  true 
undercoverage  rate  for  the  i-th  province,  and  ft  is 
treated  as  the  direct  estimator  of  0t .  Based  on  this 
transformation,  Dick  and  You  (1997)  considered  the 
following  model: 

r,=0,+e„  O^xJ/3  +  v,,  f  =  l,...,m,  (2) 
where  fi  is  a  vector  of  regression  parameters;  v;  is  a 
normal  random  effect  with  E(v,)  =  0 and  V(v,)  =  a1  „ 
where  <72is  unknown;  et  ~N(0,i//?).  The  sampling 
variance  ys?  was  treated  as  known  and  calculated  from 
a  Taylor  approximation  as 

^=(l-^)VU.)/c,2=(l-i;.)^/c,2.  The  model 
given  by  (2)  is  a  standard  application  of  Fay-Herriot 
model  for  small  area  estimation  (Fay  and  Herriot, 
1979).  However,  we  note  that  model  (2)  has  the 
following  limitations:  (i)  the  zero  mean  sampling 
model  may  not  be  true  due  to  the  nonlinear 
transformation.  The  results  obtained  from  the  model 
based  on  the  transformed  data  could  be  significantly 
different  from  those  obtained  from  the  model  based  on 
the  original  survey  estimates  (You  and  Rao,  2000);  (ii) 
the  known  sampling  variance  ysfis  a  very  strong 
assumption,  since  ysfis  indeed  a  function  of  the 
unknown  mean  0, ;  (iii)  in  model  (2)  assuming  v,.  to  be 
normal  may  not  be  appropriate  because  0t  is  the 
undercoverage  rate  which  is  between  0  and  1,  although 
similar  models  are  also  used  by  other  people,  for 
example,  Datta  et  al.  (1992)  and  Zaslavsky  (1993).  A 
more  realistic  model  may  use  log(0,)*=  xJ/3  +  vt .  To 
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oyercome  these  limitations,  we  now  propose  a  new 

model: 

Model  1: 

|»     Sampling  model:  y,.  =  m,  +  £,,   £t  -  Af(0,£2) , 

where  £f2  is  known; 
•     Linking  model: 

log(— ^— )  =  xj0  +  Vi,   Vi~N(0,a2). 

U,+Ci 

In  Model  1,  the  sampling  model  is  a  linear  model 
given  directly  by  (1)  and  is  based  on  the 
undercoverage  count.  The  sampling  variance  is 
known.  The  linking  model  is  a  log-linear  random 
effects  model  for  the  undercoverage  rate,  but  here  the 
rate  0,  is  expressed  as  a  function  of  count  m,  ,  that  is, 
0,  =  m,  /(«,  +  c, ) .  In  Model  1,  we  cannot  combine  the 
linking  model  with  the  sampling  model  to  form  a 
standard  linear  mixed  effects  model.  Thus  standard 
methods  for  linear  mixed  models  cannot  be  applied. 
This  kind  of  problems  for  small  area  estimation  is 
discussed  in  Rao  (1999)  and  You  and  Rao  (2000). 
However,  by  using  a  complete  HB  approach  with  the 
Gibbs  sampling  method,  we  can  find  the  posterior 
estimates  of  the  undercoverage  count  and  rate 
simultaneously. 

3.  GIBBS  SAMPLER 

We  now  use  Model  1  to  obtain  model-based  estimates 
of  census  undercoverage  ut .  In  a  hierarchical  Bayes 
frame,  Model  1  is  expressed  as: 


y,\u,  -NiMi4f\   i  =  l.-,m, 


(3) 


and 
u 


u,  +c 


i—\p,a2~\ogN(x]P,or2),   i  =  \,...,m,        (4) 


where  logN  denotes  a  log-normal  density  function. 
Priors  of  /Sand  a2  are  set  as:  ;r(/T)  °=  1  and 
a2  ~  IG(a,b)  with  a  and  b  are  known  positive 
constants.  Let  Y  =  (yl,...,ym)T  and  U  =(u1,.. .,um)T  , 
we  are  interested  in  the  posterior  distribution  of 
u,  given  the  data  Y,  in  particular,  posterior  mean 
E(Ui  |  Y )  and  posterior  variance  V(ut  \  Y) .  The  Gibbs 
sampling  method  (Gelfand  and  Smith,  1990)  is  used  to 
simulate  samples  for  «, .  To  implement  the  Gibbs 
sampling,  we  need  to  draw  samples  from  the  following 
full  conditional  distributions: 

.     I«||,r./»,a'j-^«P.{-^i- 


W. 


(l0g(«,  /(»,+  C,)-Xjfi)' 


2a: 


-},   i  =  l,...,m; 


•      03\Y,U,a2]~ 

m  m  •  m v 


<=1 


.=1 


M 


[a2\Y,U,P]~ 


IG{a  +  m  1 2,  b  +  £  (log ^—^-  -  xj 0) 2  1 2). 


i=i 


"'• 

Drawing  samples  from  [f3\Y,U,a  ]  and 
[C72 1 y, £/,/?]  is  straightforward.^  |  However, 
[m,  \Y,p,a2]  does  not  have  a  closed Torm.  To  draw 
samples  from  [«,  \Y,/3,a2],  Metropolis-Hastings 
updating  scheme  (see,  e.g.,  Chib  and  Greenberg,  1995) 
is  used  within  the  Gibbs  'sampler.  The  Metropolis- 
Hastings  updating  step  is  summarized  as  follows: 
Suppose  the  Markov  chain  is  at  the  k-th  iteration.   To 

update  [us\Y ,p,a2],  we  first  draw  a  candidate 
sample  m-*+1)  from  N(y,,£2),  then  with  probability 

-,/     (k)       (A+l)  %       _:_f  8\ui  'H        '"  >   n  (C\ 

a(.u\\u\  )  =  nun{  g(Mw,ft<*>,g**>)  ,1}'  ()  • 
we  accept  this  uf+V) ;  otherwise,  set  uf+1)  =  u\k) .  In 
(6),  g(-)  is  a  function  of  ut,/3,a2 given  by 


U:+C: 


g(Ui,p,cr2)  =  ^ s-exp{- 


(log(M,./(M,.+c,)-^r^)2 


2az 


■}• 


(6) 


The  estimation  of  the  posterior  distribution  of  m,  as 
well  as  £(«,.  |  n  and  V(m,  |K)can  be  based  on  the 
samples  {u\k) }  from  the  Gibbs  sampler. 

4.  IMPLEMENTATION  AND  ESTIMATION 

As  in  Dick  and  You  (1997),  we  used  1991  data  in  our 
analysis.  We  used  log-transformation  of  the  census 
count  as  the  auxiliary  variable,  that  is,  xu  =  log(c, ) , 
then         the         linking         model         for  m,  is 

log(u,  /(«,  +  c, ))  =  PQ  +  *,,/?!  +  v, .   To  implement  and 

monitor  the  convergence  of  the  Gibbs  sampler,  we 
follow  the  basic  approach  given  in  Gelman  and  Rubin 
(1992).  We  independently  simulated  L=8  sequences, 
each  of  length  t=2d,  with  d=5000.  The  first  5000 
iterations  of  each  sequence  were  deleted.  To  reduce 
the  autocorrelation  in  the  sequence,  we  took  every 
10th  iteration  for  the  left  5000  iterations,  leading  to 
500  iterations  for  each  sequence  kept  for  analysis. 
Thus  we  finally  have  L=8  sequences  with  size  n=500 
for  each  sequence. 

For  the  parameter  of  interest  ut ,  let  k^.  denote  the  j-th 
simulated  value  of  m,  in  the  1-th  sequence.  Then  the 
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m 


posterior  mean  of  m,  is  estimated  by 
ui  =  YL  X"_iM«/ / Ln.  Then  we  computed  Bl In ,  the 
variance  between  the  8  sequence  means,  w);.each 
based  on  n=500  simulated  values  of   ut ;   that  is, 

Bt/n  =  X;L=]  («,,.  -  u,  )2  /(L - 1).     Also,  let  W,  denote 

the  average  of  the  8  within-sequence  variance,  s\, 
each  based  on  n-1  degrees  of  freedom;  that  is, 
W,  =  ^  _  s?t  I L.  Then  the  posterior  variance  of  ut  is 
estimated  by  a  weighted  average  of  W,and  B, , 
namely,  &f  =  (n  -  1)W,-  /n  +  Bt/n,  where  n=500.  Note 
that  if  only  one  sequence  is  simulated,  Bi  cannot  be 
calculated.  To     estimate     the     i-th     province 

undercoverage     rate     6t ,     which     is     defined     as 

6t  =  ut  /(w,  +  c, ) ,  let  e'aj  =  u'y  /(u*,j  +  ci ) ,  then  the 
posterior  mean  of  0,can  be  estimated  by 
§i  =  X/=iSn=i^;  ^n-  Similarly,  we  can  calculate  the 
posterior  variance  of  0, .  To  monitor  the  convergence 
of  the  Gibbs  sampler,  we  calculate  VJ  =&?  +  BtILn, 

and  then  find  Rt  =V^/W,  for  each  observation.  R,  is 
known  as  a  potential  scale  reduction  factor  (Gelman 
and  Rubin,  1992).  If  fl.'s  are  near  1  for  all  of  the 
parameters  m,  of  interest,  then  this  suggests  that  the 
desired  convergence  is  achieved  in  the  Gibbs  sampler. 
In  our  study,  values  of  Rt  's  for  the  10  provinces  are 
equal  to  or  very  close  to  1,  they  strongly  suggest  that 
the  Gibbs  sampler  converged  very  well.  Table  1 
shows  the  direct  undercoverage  estimates  and  the 
posterior  undercoverage  estimates  as  well  as  the 
associated  coefficients  of  variation.  The  posterior 
estimate  of  the  census  undercoverage  count  m,  has 
smaller  coefficient  of  variation  (CV)  than  the  direct 
survey  estimate  y,  for  all  provinces  except  Ontario, 
Quebec  and  New  Brunswick.  .  Particularly  for  some 


small  provinces  such  as  P.E.I,  and  Manitoba,  the 
improvement  of  posterior  estimate  over  the  direct 
estimate  is  significant.  For  Ontario  and  Quebec,  the 
largest  two  provinces  in  Canada,  the  direct  estimate 
and  the  posterior  estimate  have  equal  CV  due  to  the 
large  sample  size  for  these  two  provinces.  An 
interesting  result  is  that  for  New  Brunswick,  the 
posterior  estimate  «,  has  larger  CV  than  the  direct 
estimate  y, ,  and  also  «,  is  quite  different  from  y,. . 
This  is  more  obvious  in  terms  of  undercoverage  rate, 
where  directly  estimated  rate    rt   is  3.25%  and  the 

posterior  rate  0,  is  only  2.55%  for  New  Brunswick. 
Actually  province  New  Brunswick  is  an  outlier  in  the 
model.  More  detailed  analysis  will  be  given  in  the 
next  section. 

5.  TEST  OF  MODEL  FITTING 

To  test  the  overall  fit  of  Model  2,  we  use  the  method 
of  posterior  predictive  p  value  (Meng,  1994;  Gelman 
et  al.,  1995).    Let  T(y,0)  be  a  discrepancy  measure 

depending  on  the  data  y  and  on  parameters  6 .  Let  9*  v 
represent  a  draw  from  the  posterior  distribution  of  8 : 

and  let  y*  represent  a  draw  from  f{y\9*),  then 
marginally  y*  is  a  sample  from  the  posterior 
predictive  distribution  /(y|yofa),  where  yobs 
represents  the  observed  data.  The  posterior  predictive 
p  value  is  defined  as 

p  =  Pr(ny',e)>nyobs,d)\yobs).  Note  that  the 
probability  is  with  respect  to  the  posterior  distribution 
given  the  observed  data. 

This  is  a  natural  extension  of  the  usual  p  value  in  a 
Bayesian  context.  If  a  model  fits  the  observed  data, 
then  the  two  values  of  the  discrepancy  measure  are 
similar.  In  other  words,  if  the  given  model  adequately 
fits  the  observed  model,  then  T(yobs,0)    should  be 


Table  1.  Census  undercoverage  estimation 


Province ' 

yt    cvty,) 

u,      CV(h,.) 

rt% 

CV(f,) 

e;% 

cvw,) 

NFLD  , 

11566   0.16 

10782   0.14 

1.99 

0.16 

1.86 

0.13 

PEI  '    ; 

1220  0.30 

1486   0.19 

0.93 

0.30 

1.13 

0.19 

NS 

17329   0.20 

17412   0.14 

1.89 

0.20 

1.90 

0.14 

NB 

24280   0.14 

18948    0.17 

3.25 

0.13 

2.55 

0.17 

QUE      ;■ 

184473   0.08 

189599   0.08 

2.58 

0.08 

2.65 

0.08 

ONT  ,! 

381104,  0.08 

368424   0.08 

3.64 

0.08 

3.52 

0.08 

MAN    ■ 

,■.••20691    0.2-1- 

21504   0.14 

1.86 

0.20 

1.93 

0.14 

SASK 

18106   0.19 

18822   0.14 

1.80 

0.18 
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near  the  central  part  of  the  histogram  of  the  T(y*,0) 

values  if  y*  'is  generated  repeatedly  from  the  posterior 
predictive  distribution.  Consequently,  the  posterior 
predictive  p  value  is  expected  to  be  near  0.5  if  the 
model  adequately  fits  the  data.  Extreme  p  values  (near 
0  or  1)  suggest  poor  fit.  Computing  the  p  value  is 
relatively  easy  using  the  posterior  simulation  from  the 
Gibbs  sampler.  For  each  simulated  value  0* ,  we  can 
simulate  y*  from  the  model  and  compute  T(y' ,0*) 
and  T(yobs,0*) .  Then  the p  value  is  estimated  by  the 
proportion     of     times     that      T(y*,0*)      exceeds 

T(yobs,0*).      In    the    present   context    of   census 

undercoverage,  the  discrepancy  measure  that  we  used 

for  overall  fit  is  T(y,u)  =  ^l.(yi  -ut)2 /£?  .    Similar 

discrepancy  measure  is  also  used,  for  example,  in 
Datta  et  al.  (1999).  We  computed  the  p  values  for  all 
8  parallel  runs.  Under  Model  1,  the  smallest  p  value  is 
0.374,  the  largest  one  is  0.482,  and  the  average  p  value 
is  0.433.  Thus  we  have  no  indication  of  any  lack  of 
overall  fit.  The  p  value  strongly  suggests  the  adequacy 
of  the  model. 

We  also  computed  two  statistics  that  are  useful  in 
assessing  model  fit  at  the  individual  observation  level. 

First,  we  computed  p'  =  Pr(y*  <  yi{obs)  |  yobs) ,  which 

provides  information  on  the  degree  of  consistent 
overestimation  or  underestimation  of  the  observed 
data.  Second,  we  computed 

r*=E*( yt  -  yi(obs)  |  yobs ) / yJVar(yi  \yobs)  ,  the 

standardized  predictive  residual.  So  we  used  the  mean 
of  the  predictive  distribution  for  yt  as  the  center  and 

the  variance  as  a  measure  of  the  spread.  These 
standardized  residuals  can  be  used  to  compare  model 
fit  for  individual  observations.   We  used  the  posterior 

simulations  y*  to  compute  p*  and  r* .  Under  Model 

2,  p*  ranges  from  0.28  to  0.87,  with  a  mean  of  0.51 
and  a  median  of  0.49,  indicating  no  consistent 
overestimation  or  underestimation.     Among  the   10 

provinces,  PEI  has  p*  =0.28  and  New  Brunswick  has 

p*=0.87,  indicating  that  the  model  may  not  be 
appropriate , for 'these  two  provinces,  especially  for 
New  Brunswick.  For  the  standardized  residuals  r* , 
New  Brunswick  has  the  largest  positive  a;* =1.13, 
which  indicates  that  the  current  model  underestimates 
the  census  undercoverage  count  with  respect  to  the 

observed  value.  PEI  has  the  largest  negative  r*  =-0.57, 
which  indicates  that  the  model  overestimates  the 
undercount  compared  to  the  observed  value.   In  terms 


of  absolute  value  of  /;.*-,  New  Brunswick  has  the 
largest  standardized  residual,  p'and  r*  can  Also  'be 
used  as  criteria  for  outlier , detection.  An'  outlier  is 
defined  as  an  observation  "that  is  extreme  (in  the  value 
of  some  statistic)  relative  to  its  predictive  distribution 
under  the  model.  Because  outliers  represent  some  data 
patterns  that  are  not  well  described  by  the  model,  they 
suggest  areas  in  which  modifications  of  the  model 

should  be  considered.  In  terms  of  the  values  of  p*  and 
r* ,  New  Brunswick  is  obviously  an  outlier.  Dick  and 
You  (1997)  also  obtained  similar  findings.  Since  New 
Brunswick  always  appears  to  be  an  outlier  in  various 
models,  we  suggest  that  further  work  should  be  done 
for  the  direct  survey  estimate  of  the  census 
undercoverage  for  New  Brunswick. 

6.  CONCLUSION  REMARKS 

In  this  paper,  we  have  proposed  a  nonlinear  mixed 
model  for  the  problem  of  census  undercoverage 
estimation.  The  model  involves  a  linear  sampling 
model  and  a  log-linear  linking  model.  We  cannot 
combine  the  sampling' model  with  the  linking  model  to 
form  a  usual  linear  mixed  model  to  which  standard 
methods  can  be  applied.  However,  with  the  Gibbs 
sampling  method  and  Metropolis-Hastings  algorithm, 
we  find  the  HB  estimator  of  the  census  undercoverage 
numerically.  Our  data  analysis  shows  that  the 
proposed  model  fit  the  data  quite  well.  The  HB 
method  improves  the  direct  survey  estimates  and  also 
can  be  used  to  identify  possible  outliers  in  the  data. 
This  paper  itself  also  provides  a  good  example  for  the 
case  of  unmatched  sampling  and  linking  models  in 
small  area  estimation  as  discussed  in  You  and  Rao 
(2000). 
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